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DISPLACEMENT MODELS FOR THUNDER ACTUATORS HAVING GENERAL LOADS 

AND BOUNDARY CONDITIONS 

ROBERT WIEMAN*, RALPH C. SMITH*, TYSON KACKLEY*, ZOUBEIDA OUNAIES§, AND JEFF BERND*! 

Abstract. This paper summarizes techniques for quantifying the displacements generated in THUN- 
DER actuators in response to applied voltages for a variety of boundary conditions and exogenous loads. 
The PDE models for the actuators are constructed in two steps. In the first, previously developed theory 
quantifying thermal and electrostatic strains is employed to model the actuator shapes which result from the 
manufacturing process and subsequent repoling. Newtonian principles are then employed to develop PDE 
models which quantify displacements in the actuator due to voltage inputs to the piezoceramic patch. For 
this analysis, drive levels are assumed to be moderate so that linear piezoelectric relations can be employed. 
Finite element methods for discretizing the models are developed and the performance of the discretized 
models are illustrated through comparison with experimental data. 

Key words. THUNDER actuators, shape model, displacement model 

Subject classification. Structures and Materials 

1. Introduction. THUNDER actuators offer the capability for generating large strains and forces due 
to a variety of mechanisms including improved robustness through the manufacturing process and increased 
electromechanical coupling due to their inherent shape. However, the full capabilities of these actuators 
have not yet been completely quantified either experimentally or analytically due to their relatively recent 
genesis and the fact that their behavior differs quite substantially from standard unimorphs or bimorphs. 
In this paper, we discuss modeling techniques for quantifying the displacements generated by THUNDER 
actuators in response to applied voltages for a variety of boundary conditions and exogenous loads. The 
development of corresponding finite element techniques is also addressed and the accuracy of the resulting 
finite dimensional models is illustrated through comparison with experimental data. 

Model development is considered in two steps: (i) the characterization of the actuator shape as a function 
of the manufacturing process and (ii) the development of a PDE model for the actuator behavior based on 
Newtonian principles. The first component has been addressed in previous investigations [3] and only those 
details necessary for the development of the subsequent PDE model will be discussed. As detailed in [3], 
the characteristic curved shape of THUNDER actuators is due primarily to differing thermal coefficients 
in the constituent materials, which produce thermal stresses in the combined actuator during cooling, and 
secondarily to the reorientation of dipoles during repoling. We note that the quantification of strains due to 
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thermal gradients has been investigated for a variety of applications (e.g., see [2, 4, 5, 6]) with certain aspects 
having been considered for THUNDER actuators [3, 7]. The quantification of stresses and strains due to 
repoling is based on domain theory for general ferroelectric materials [8, 9, 12, 13]. Thin shell theory is then 
employed to develop PDE models which quantify the stresses and displacements throughout the actuator 
when voltage is applied to the piezoceramic patch. For this analysis, it is assumed that the actuators are 
operating at low to moderate drive levels for which linear piezoelectric relations are adequate. Techniques 
for extending these models to regimes in which the piezoelectric response is nonlinear and hysteretic are 
under investigation and will utilize methods outlined in the concluding remarks. 

Because the PDE model is infinite dimensional, approximation techniques must be considered to obtain 
a finite dimensional model which is appropriate for implementation. This is accomplished through a hybrid 
finite element approach utilizing linear and cubic Hermite basis functions. This numerical approach differs 
from that employed in [14], where a NASTRAN model was employed to predict dome heights, in that the finite 
element method was developed directly for the PDE model used to characterize the physical mechanisms for 
the actuator. Hence the resulting finite dimensional system incorporates the physical properties associated 
with the differing constituent materials thus permitting a detailed analysis of various aspects of the actuator 
dynamics (e.g., stresses or strains at various points along the length of the actuator). This approach also 
permits direct extension of the numerical method to nonlinear structural models for high drive level dynamics 
as well as models which incorporate the hysteresis and constitutive nonlinearities inherent to piezoceramic 
materials at moderate to high drive levels. 

The manufacturing conditions for THUNDER are outlined in Section 2 and a model which quantifies 
the resulting curved shape is summarized in Section 3. The PDE model quantifying the displacements is 
then presented in Section 4 along with boundary conditions which characterize a variety of experimental 
setups. The numerical approximation techniques are discussed in Section 5 and examples illustrating the 
performance of the resulting finite dimensional model are presented in Section 6. Finally, future work, 
including techniques to extend the model to nonlinear and hysteretic regimes will be outlined in Section 7. 

2. Actuator Geometry. THUNDER actuators are typically comprised of a piezoceramic wafer, a 
metallic backing material, hot melt adhesive layers, and an optional metallic top layer as depicted in Fig- 
ure 2.1a. As detailed in [3], materials commonly employed for backing layers include aluminum, stainless 
steel and brass while LaRC-SI is employed as the adhesive. 

During the manufacturing process, the materials are placed in a vacuum bag and heated to 325° C under 
a pressure of 241.3 kPa. During the cooling process, the LaRC-SI solidifies at approximately 270° C and 
subsequent cooling produces curvature in the actuator due to differing thermal coefficients of the constituent 
materials. Because the Curie temperature for PZT-5A (350° C) is in the proximity of the manufacturing 


LaRC-SI Film 


LaRC-SI Film 



Fig. 2.1. (a) Components of a THUNDER actuator; (b) Curvature observed in a THUNDER actuator. 
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temperature, the final step in the fabrication process is comprised of repoling the material through the 
application of a sustained DC voltage. 

As illustrated in Figure 2.1b, THUNDER actuators have a characteristic dome shape due to the man- 
ifestation of differing thermal properties in the PZT and backing material during the cooling process and 
the rotation of dipoles during repoling. In general, curvature will occur in both component directions in a 
rectangular actuator; however, for the models developed here, we consider actuators whose width is small 
as compared with the length so that motion is predominantly in one dimension. Finally, we note that the 
curvature in the actuators occurs only in regions covered by the piezoceramic patch and the end tabs remain 
straight. 

3. Model for the Manufacturing Process. A necessary step before developing a model which 
predicts the actuator displacements under various drive conditions is the quantification of the curvature 
produced in the manufacturing process. We summarize here the characterization of the stresses produced 
during cooling and repoling which in turn produce the curvature. Details regarding this component of the 
model can be found in [3]. 

To accommodate various constructions, we consider actuators with N layers and consider a coordinate 
system in which the x — z plane corresponds with the outer edge of the backing material and the ^-coordinate 
extends through the thickness of the actuator. The width of the jtli layer is denoted by bj while hj indicates 
the thickness of each layer as depicted in Figure 3.1. The Young’s modulus and thermal coefficient for the 
jth layer are respectively denoted by Ej and aj. The strain at the outer edge of the backing material ( y = 0) 
is denoted by £o, and k denotes the curvature at the neutral axis. The change in temperature during the 
bonding process is indicated by AT. To incorporate the strains due to repoling, it is also necessary to employ 
the Poisson ratio v and saturation electrostriction A s for PZT-5A. 

As detailed in [3], the balancing of forces and moments due to thermal and electrostatic stresses yields 
the linear system 


(3.1) 


A£ = f 


where £ — [eq, k] t and 


(3.2) 


A = 


E;=i Ejbj (hj - hj _! ) -A Ej=i Ejbj (hj - h)_ x ) 

\ Ejir Ejbj(h) - h)_ x ) £f=i Ejbj(h) - hj_ x ) \ 


f = 


£?= i Ejbj(aj AT - 3/26v\ 8 )(hj - hj-,) 

\ Ej=i Ejbj (a j AT - 3/2SuX s )(hj - hj_ x ) 


The Kronecker delta, defined by 



1 , if y is in the piezoceramic layer 
0 , otherwise 


isolates the electrostatic strains due to repoling to the piezoceramic layer. 

To solve for £q and k, and hence obtain the final radius of curvature R = 1/k, it is necessary to 
obtain values for Ej,otj for each of the constituent layers in addition to determining v\ s for the PZT 
compound being employed. While the Young’s modulus and thermal coefficients are catalogued for various 
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Fig. 3.1. Orientation of the composite THUNDER actuator with five layers. 


PZT compounds and various metallic backing materials, they are temperature-dependent, and not easily 
quantified, for the LaRC-SI. Hence these parameters are typically estimated through a least squares fit 
to data for the constituent materials and specific manufacturing conditions under consideration. Details 
regarding the validity of the model for a variety of materials and geometries can be found in [3]. 

4. Displacement Model. The system (3.1) quantifies the radius of curvature R = 1/k for differing 
constituent materials, material dimensions, and manufacturing conditions. In this section, we develop models 
quantifying displacements produced in the actuators through the input of voltages to the PZT or applied 
loads. Because end conditions crucially affect the measured displacements, we consider a variety of boundary 
conditions. Finally, we consider low drive regimes in which the relation between applied voltages and 
generated strains are approximately linear with minimal hysteresis so that linear piezoelectric relations can 
be employed. 

When modeling the actuator, we consider tw r o regimes. In the first, the entire actuator (including 
the tabs) is assumed to have the same initial curvature. This approximation to the geometry significantly 
simplifies the numerical implementation of the model but imposes the assumption that the tabs are initially 
curved. In the second configuration, only that region covered by PZT is assumed to be initially curved and 
the end tabs are taken to be flat in the absence of an applied voltage or load. This accurately represents the 
initial configuration of the actuator after the manufacturing process. One of the objectives when validating 
the model is to compare the performance of both models and ascertain tab dimensions when the latter model 
is sufficiently accurate. 

4.1. Model 1. We consider initially the model which results from the assumption that the tabs have 
the same curvature as that portion of the actuator covered by the PZT patch. The radius of curvature is 
denoted by R (recall that R = 1/k can be predicted using the model summarized in Section 3) and the 
backing material is assumed to have width 6, thickness h, and Young’s modulus E. The corresponding 
parameters for the PZT layer are respectively denoted by b pe ,h pe and E pe . The backing material is assumed 
to extend from 9 = 0 to 9 =' L and the region [71 , j 2 ] covered by the patch is delineated by the characteristic 
function 

, . f 1 , 7i < 7 < 72 

(4-1) Xpe = < 

( 0 , otherwise 

where 71 = R9\ and 72 = R9 2 and 61 , d 2 denote the angles subtended by the patch. The longitudinal and 
transverse displacements, which are coupled due to the curvature, are respectively denoted by v and w. 

Under the assumptions of linear displacements, negligible rotational effects and shear deformation, and 
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linear stress-strain relations, force and moment balancing yields the static equations 

1 dNg _ _ 

(4.2) R d9 _ 

1 <fM e 1 . 1 (PM e 

R 2 d6 2 + R Ae ~ Qn R 2 dd 2 ' 


Here Ng and Mg denote the internal force and moment resultants and Mg represents the external moment 
generated by applied voltages to the patches. Finally, qg and q n respectively denote applied longitudinal 
and normal loads to the actuator. As detailed in [1], the modeling equations (4.2) are consistent with the 
restriction of modified Donnell-Mushtari shell equations to the actuator geometry. 

The internal resultants incorporate the material properties of the backing material and PZT and, as 
derived in [1], are given by 


(4.3) 


„ y , 1 dv 

N ’ = Ehl .Rde + 


' dv 
MB 


, , -Eh 3 cPw „ 

' 0 ~ 12 R 2 dd 2 + pe 


i)+M£( 

a.2 f dv \ 
2 R \d0 +W ) 


+ w — 


02 <Pw 
2 R 2 d6 2 


Xpe(9) 


<*3 


drw 1 


3 R 2 d6 2 


X P e(0) 


where a 2 = (h/2 + h pe ) 2 - (/i/2) 2 and 03 = (h/2 + h pe ) 3 - (h/2) 3 . 

For low to moderate drive regimes in which the linear piezoelectric equations are sufficiently accurate, 
the external moment generated by the PZT in response to an applied voltage V is 


(4.4) 


Mg = Ep ' d31 V(h + h pe ) X pe{6) 


where d^i is the linear piezoelectric constant. 

We note that the model (4.2), with resultants given by (4.3) and (4.4), neglects material contributions 
due to the LaRC-SI. If desired, these contributions can be incorporated in the manner described in [1], We 
also note that in the strong form (4.2), differentiation of the discontinuous material parameters and patch 
inputs yields unbounded components in the model. This necessitates the consideration of an appropriate 
weak form of the model. As a prelude, however, it is necessary to specify appropriate boundary conditions. 

We consider four sets of boundary conditions which model the constraints commonly employed in ex- 
periments: fixed-end, pinned-end, sliding-end and free-end conditions. These boundary conditions can be 
applied at either end of the beam; to simplify the discussion, we summarize them at the left end (9 = 0) and 
note that similar expression hold at 9 = L. 


(i) Fixed-End Conditions 


(4.5) 


w(0) = 0 

-( 0 ) = >)=0 


(ii) Pinned-End Conditions 


(4.6) 


w(0) = 0 

w(0) = Mg{ 0) = 0 
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(iii) Sliding-End Conditions 


w( 0) = u(0) tan(<^>i) 

(4.7) M e { 0) = 0 

Ne( 0) = -Qe(O) tan(</> c ) 


(iv) Free-End Conditions 


(4.8) Ar 0 (O)=M„(O)=O 

In the sliding end condition, Qg denotes the shear force resultant and (pi, <p c respectively denote the 
initial angle of the actuator and the angle obtained after a load is applied (see Figure 4.1). It can be noted 
that to within a first-order approximation, (pi and <p c are related by the expression <p c — <Pi + ■ Finally, 

for implementation purposes, it has been observed that physically reasonable results can be obtained with 
the approximation Qg = 0 which is typically enforced in first-order shell theory. This yields the natural 
boundary condition Ng = 0 which is easily implemented. 

We note that care must be exhibited when specifying boundary conditions at the left and right ends of 
the actuator to ensure model well-posedness. For example, the specification of free-end conditions at both 
x — 0 and x — L will yield rigid body modes and hence will not enforce unique solutions since solutions 
differing by a constant will be equivalent. For the experiments reported in Section 6, fixed-end conditions 
were enforced at x = 0 and sliding-end conditions were employed at x = L. 

To accommodate the discontinuities due to the patch and to reduce smoothness requirements on the 
basis functions employed for numerical approximation, we consider corresponding weak forms of the modeling 
system. The state space is taken to be X = L 2 (U) x L 2 (Q) where fi = [0, L\. The test functions depend 
upon the boundary conditions under consideration. For fixed, pinned, or sliding-end boundary conditions at 
0 = 0 and free-end conditions at 6 = L, we respectively employ the spaces 

v = {(<P,<f) € H l x H 2 I #)) = 0,^(0) = v>'(0) = 0} 

(4.9) V = {(0 >¥ ») € H 1 x H 2 \ cP( 0) = 0,^(0) = 0} 

V = {(<P,<p) e H 1 x H 2 | tp( 0) = <£(0)tan((fo)} . 

We note that the constraints Mg(L) = N$(L) = 0 are natural boundary conditions which do not require 
any restriction of the underlying Sobolev spaces. Analogous definitions are employed when considering other 
combinations of boundary conditions. 

A weak form of the model is then 

f 0 

(4 - 10) 

f { W ~ ™ + m = 0 

for all {(p, <p) in the appropriate space V. 
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4.2. Model 2. A second model is derived under the assumption that the actuator region covered by 
the patch is curved while the tabs are initially straight. This model is constructed by coupling flat and 
curved beam models through appropriate interface conditions. 

To specify the geometry, we let 7 denote the arclength with 7 = 0 at the left end of the actuator 
as depicted in Figure 4.1. We assume that the tabs have equal length t and that the portion covered by 
the piezoceramic patch has arclength s = (62 — 0\)R. The region covered by the patch is denoted by 
fl pf , = [71,72] = [f,s + <] while fi = [0, L] again denotes the domain for the full actuator. The characteristic 
function Xpe , defined in (4.1) delineates that portion of the structure covered by the patch. For the tabs, 
the arclength is designated by 7 = x whereas it has the form 7 = R0 in the curved portion of the actuator. 

For the curved portion of the actuator, force and moment balancing again yield the coupled relations 


(4.11) 


1 dNg _ . 

R dd 9e 

1 (PMg 1 „ 1 (PMg 

R 2 dff 1 + R 1(l ~ <ln R ' 2 dff 2 ' 


where Ng,Mg and Mg are defined in (4.3) and 4.4). The tabs have infinite radius of curvature which yields 
the uncoupled relations 


(4.12) 


dN x _ „ (PM X _ . 

dx 9x ’ dx 2 9n 


where the resultants are given by 

dn 

(4.13) N x = Eh— 


M x 


Eh 3 cPw 
~l2d^' 


Finally, the displacements and slopes at 71 and 72 are coupled through the interface constraints 


(4.14) 


lim v(j) - lim 11(7) 
7~>7f 7->7! + 


lim u(7) = lim ^(7) 

7->7^ 7-+7^ 


lim w( 7) = lim w( 7) 
7->7f 


lim w(-f) = lim 10(7) 

7->72 7 -> 7 ? 


lim _ 7^(7) = lim ^(7) 

7— dx 7—1 ►7^’ dx 


.. dw , x .. dw 
lim —(7) = lim — 
7 — >7” dx — ^7+ dx 


(j)- 


Formulation of the model in terms of the arclength 7 yields the unified relations 


(4.15) 


dN-y _ , 


d 2 M-y 1 j . (PM-, 

-^f L+ R N ^ qn ~^f L - 
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with 


.. f dv w\ dv h pe w a 2 (Pw 1 , 

7 (^7 + RJ + pe [ hpe d^ + "IT “ T Xpeh) 


—Eh 3 drw la^dv a^w 0,3 d 2 wl . 

M -< = —i? + E >- It ^ + 2K - y w\ Xreh) 


= l -V(h + h rt ) Xp ,h)- 

The radius of curvature is taken to be R = 00 for the tabs. When employing the formulation (4.15), it 
must be noted that the second derivatives may not exist at the points 71 and 72. Finally, fixed, pinned and 
sliding-end boundary conditions are enforced by employing the constraints (4.5) - (4.7) at 7 = 0 or 7 = L. 

The construction of the weak model formulation is analogous to Model 1 with the exception that second 
derivatives of the transverse test functions may not exist at the interface points [71,72]- To illustrate, the 
space of test functions for an actuator with fixed-end conditions at 7 = 0 and sliding-end conditions at 7 = L 
is 


V = | (<j>,<p) e H 1 x H 2 I 0(0) = 0,¥>(0) = v>'(0) = 0 ,<p(L) = 0(L)tan(0<) and 

lim ip"( 7) 7^ lim ip"( 7), lim </?"( 7) 7^ lim <^"(7) 1 . 

7-nr 7-nJ ■y-'tt ) 

Analogous spaces are employed for the remaining combinations of boundary conditions. A weak form of the 
model is then 

(4.17) 

for all (0, ip) £ V. 

5. Numerical Approximation Techniques. To approximate the solution of (4.10) or (4.17), we 
consider Galerkin techniques with basis functions chosen to satisfy smoothness requirements as well as 
boundary and interface conditions. We consider first the system which arises when discretizing the model 
for the uniformly curved actuator. 


5.1. Model 1. We consider Galerkin approximation for v and w which are respectively based on linear 
and cubic Hermite functions. To define the bases, consider a uniform partition of [0, L] with gridpoints 
9 i = ill, h = L/N, i — 0, • • • , N. For i = 1, • • • , N — 1, linear splines are taken to be 

(0-0<-i) , oe^A] 

(5-i) 0i(0) = ^< (e i+1 -0) , ee[0 i ,e i+ 1 ] 

0 , otherwise . 

The cubic Hermite basis functions employed to specify w and w' are given by 


(0-O i - 1 ) 2 [2(9 i -6) + h] , 0G[0i_i ,0i] 

^0(6) = ^ (6 i+1 -0) 2 [2(0 i+1 -0)-h] , 9 £( 0 i,e i+ 1] 

0 , otherwise 
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and 


(0 - 0i-i) 2 (e - 0i) , 

<pa(0) = \ (o i+ i-ey 2 (e-0i) , 9e(0i,9 i+1 ] 

0 , otherwise . 

As detailed in [11], all three sets of basis functions vanish identically outside the interval [0 4 -_i,# j+ j]. The 
basis functions </>o,<A)o><A)i and 4>n,ipno,<Pni are defined similarly on the intervals [0, 9\\ and [0 jv_i,0jv] 
(see [11, pages 49 and 57]). The essential boundary conditions (i)-(iii) summarized in (4.5)-(4.7) are enforced 
by omitting or forming linear combinations of the boundary basis functions. 

The displacements are then represented as linear combinations of basis functions with coefficients de- 
termined by enforcing the constraints provided by the weak form of the model. To illustrate, consider the 
discretization of the model satisfying fixed-end conditions at 9 = 0 and pinned-end conditions at 9 = L. The 
approximate solutions are taken to be 

N - 1 

V N {0) = Y V 3 ( t > ^ e ) 

i=i 


w n { 9) = Y Wjtpj o(0) + Ytijtpji(9) 
l=i l=i 

in the subspace H N C V. Analogous expressions are employed for other combinations of boundary condi- 
tions. 

A matrix system is obtained by considering the approximate solution in (4.10) with basis functions 
employed as test functions (this is equivalent to projecting the system (4.10) onto the finite dimensional 
subspace H N ). This yields the linear relation 

(5.3) Kv = f 

where v = [iq , • ■ • , w/v_i ,w i , - • • ,wn - i , , ■ ■ ■ , wn] denotes the vector of unknown coefficients. 

5.2. Model 2. The formulation of approximation techniques for Model 2 is accomplished in a similar 
manner. In this case, we consider uniform partitions on each of the subintervals [0, 71], [71,72], [72, £] arid 
define basis functions on each subdomain. The displacements v N = [u 1 /v , Vj] and w N = f«>[ v , w^] 
are defined in a manner analogous to (5.2) with the interface conditions enforced through the constraints 


«i V (7i) = 1^(71) 


» v 2 (72) = (72) 


wf f (7i) = W'f(7i) . ^(72) = Ws( 72) 


dw¥ . dwn . , dw n 

• tr 


dwjs . 

(72) = -37- ( 72 ) 


Orthogonalization against the test functions then yields a corresponding linear system 

Kv = f 


which can be solved to obtain the displacement and slope coefficients at each of the nodes in the subdomains. 
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6. Experimental Validation. In this section, we illustrate the performance of the combined model 
through comparison with experimental data. In the first example, the strain model outlined in Section 3 is 
used to predict the shape of a variety of actuators as a function of the manufacturing process. This provides 
the initial geometry employed in Section 4 when characterizing displacements generated through voltage 
inputs to the patches. The performance of the displacement model is illustrated in Example 2. 

6.1. Example 1. Actuator Shape 

The model summarized in Section 3 quantifies the thermal and electrostatic strains and resulting changes 
in curvature generated during the cooling and repoling of the actuator during the manufacturing process. 
This provides a means of characterizing the radius of curvature and dome height for an actuator as functions 
of properties of the constituent materials as well as the dimensions of these materials. We consider in this 
example an actuator construction comprised of a stainless steel bottom layer, LaRC-SI, PZT-5A and a 
protective LaRC-SI top layer (there is no metallic top layer). The width of all materials was 0.5 inches and 
the PZT was 1.5 inches in length while the stainless steel was 2.5 inches in length. Hence s = 1.5 and t = 0.5 
in Figure 4.1. The PZT was 8 mils thick while the LaRC-SI had a mean thickness of 1 mil. Actuators were 
constructed with steel thicknesses ranging from 1 mil to 20 mils to illustrate the effect of backing material 
thickness on the final dome height of the actuator. Note that the dimensions of the THUNDER devices 
considered here permit the use of the 1-D model. 

The dome heights h predicted by the relation (3.1) are compared with experimental data from actuators 
having steel thicknesses ranging from 3 mils to 20 mils in Figure 6.1. In the model, the parameter values 
E pe = 177 x 10 9 N/m 2 , E S i = 7.45 x 10 9 N/m 2 , E stee i = 173 x 10 9 N/m 2 , a pe = 0.8 x HT 7 , a S i = 46 x HT 6 
(23-150° C), asj = 60 x 10~ 6 (23-150° C), a stee i = 9.8 x 10 -6 , A s = 0.6 x 10 -3 and v = 0.3 were estimated 
through a least squares fit to the full data set. It is observed that the model predicts both the trends and 
magnitudes for the dome heights to within 5% relative accuracy for all steel thicknesses except 6 mils for 
which the relative error was 8%. This provides a characterization of the actuator shape which can then be 
employed when modeling subsequent displacements or forces. 
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Fig. 6.1. Model predictions and experimentally measured dome heights h as a function of steel thickness. 
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6.2. Example 2:. Actuator Displacements 

To illustrate the performance of the displacement model presented in Section 4, we consider the con- 
struction described in Example 1 for actuators having fixed-end conditions at the left edge (7 = 0) and 
sliding-end conditions at the right edge (7 — L). Model 2 was used to predict the displacements generated 
at the patch center (7 = t + s/2) for a variety of steel thicknesses and input voltages. 

For fixed voltage levels, model predictions for the displacement to voltage ratio as a function of steel 
thickness are compared with experimental data in Figure 6.2a. It is observed that the model accurately 
predicts the displacement for actuators having steel thicknesses of 3 mils, 8 mils and 10 mils with discrepancies 
observed at 1 mil and 6 mils. We note that the LaRC-SI has the same thickness as the steel at 1 mil and 
we hypothesize that in this case, unmodeled viscoelastic properties of the LaRC-SI may be dominating the 
elastic properties of the steel. The error in the 6 mil prediction reflects the discrepancy observed in the dome 
height prediction for that thickness. 

To further illustrate the performance of the model at low drive levels, the predicted displacements for 
input voltage levels of 20 V, 80 V and 120 V are compared with experimentally measured displacements in 
Figure 6.2b. The steel thickness in this case was 10 mil. It is observed that within this (approximately) 
linear range, the model accurately predicts the displacement for a variety of input levels. 

7. Concluding Remarks. The model described here provides a technique for quantifying both the 
initial shape of THUNDER devices due to the manufacturing process and displacements generated by ap- 
plied voltages. The actuator shapes were modeled through the quantification of thermal and electrostatic 
strains while Newtonian principles were used to derive PDE models characterizing displacements for a va- 
riety of boundary conditions and exogenous loads. Both components of the model were illustrated through 
comparison with experimental data. 

It should be noted that linear theory was employed when deriving both components of the model and 
degradation of performance is expected in high drive regimes. The extension of both components to nonlinear 
regimes is under current investigation. The extension of the thermal model to include nonlinear effects is 
being considered in the context of theory in [10]. To incorporate the constitutive nonlinearities and hysteresis 
inherent to piezoceramic materials at moderate to high drive levels, the models developed in [12, 13] are 
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Fig. 6.2. (a) Model predictions and measured displacements as a function of steel thickness; (b) Modeled and measured 
displacements as a function of input voltage: xxxx (Data), (Model). 
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being combined with the current model to accommodate large inputs. Finally, nonlinear shell theory will be 
employed to ascertain limitations in the linear PDE presented here when displacements are large. 
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